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Abstract 

Atomic scale computer simulations using density functional theory were used to investi¬ 
gate the behaviour of tin in the tetragonal phase oxide layer on Zr-based alloys. The Sn^. r 

II 

defect was shown to be dominant across most oxygen partial pressures, with S n Zr charge 
compensated by Vq* occurring at partial pressures below 1CP 31 atm. Insertion of addi¬ 
tional positive charge into the system was shown to significantly increase the critical partial 
pressure at which S n Zr is stable. Recently developed low-Sn nuclear fuel cladding alloys 
have demonstrated an improved corrosion resistance and a delayed transition compared to 
Sn-containing alloys such as Zircaloy-4. The interaction between the positive charge and 
the tin defect is discussed in the context of alloying additions such as niobium and their 
influence on corrosion of cladding alloys. 


1 Introduction 

1.1 Corrosion of zirconium 

Zirconium alloys are used in pressurised water reactors due to their low capture cross-section 
for thermal neutrons, high corrosion resistance and acceptable mechanical properties. Initially, 
corrosion follows an approximately cubic rate law [Tj . After a few microns of oxide growth, there 
is a sudden increase in the corrosion rate as the protective layer breaks down; this process is 
generally referred to as transition, or break-away corrosion. Post-transition, a reduced corrosion 
rate is again observed. This process repeats in a cycle; the time from initial corrosion to transition 
and subsequent later transitions is also composition dependent m- 


i 


During corrosion hydrogen is produced, which arrives at the metal oxide interface. Hydrogen 
has a low solubility in Zr 02 so that after recombining with electrons, the hydrogen atoms move 
either into the coolant, or into the oxygen-saturated a-Zr beneath the oxide (in which hydrogen 
exhibits a much higher solid solubility [3]). When the solubility limit for hydrogen in a-Zr is 
exceeded, zirconium hydrides precipitate. The fraction of hydrogen produced during corrosion, 
that moves into the cladding metal, is called the hydrogen pick-up fraction (HPUF) and is alloy 
dependent. It has been suggested [3] that alloys which exhibit a lower HPUF have a more 
electrically conductive oxide layer, allowing ion-electron recombination to occur further from 
the metal-oxide interface thereby lowering the probability of the hydrogen entering the cladding 
metal. A lower HPUF is desirable, since there is a regulatory upper limit to the amount of 
hydrogen allowed in the cladding metal [5] . 

1.2 Oxide layer composition 

Various experimental techniques have demonstrated that the pre-transition oxide layer formed on 
the surface of zirconium alloys is composed of several distinct regions. From the outside inwards: 
an outer layer of predominantly columnar monoclinic Z 1 O 2 , a 50-80 nm layer containing equiaxed 
tetragonal ZrC> 2 , often a 100-200 nm layer of sub-stoichiometric Zr-oxide and finally a-Zr metal 
(saturated with oxygen close to the nretal-oxide interface) [6H8]. 

Early studies of the sub-stoichiometric oxide layer were unable to identify the crystal structure 
other than it resembled distorted a-Zr [9]. Subsequent TEM work hinted at several possible sub- 
stoichiometric zirconia structures, such as cubic ZrO [10][TTj or an ordered oxygen solid solution 
within w-Zr with an approximate stoichiometry of Z^O [121 . Random structure searching of 
the Zr-system using density functional theory have identified previously unknown energetically 
stable structures USE!, recent studies have identified a metastable ZrO structure, with an HCP 
lattice similar to a-Zr Hsus. 

It is generally agreed that an increased proportion of tetragonal Zr02 phase is found adjacent 
to the metal-oxide interface, where it is stabilised by a combination of grain size, compressive 
stress and the presence of alloying (i.e. dopant) elements [T7]-Q2!]- Immediately post-transition, 
oxides have been shown to contain very little tetragonal or sub-stoichiometric zirconia away 
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from the metal oxide interface, suggesting that these phases are transformed to monoclinic ZrC >2 
during transition [5], 

At the metal oxide interface, the oxygen partial pressure is expected to be extremely low, 
however the limited amount of atom probe [20] and TEM m studies thus far have demonstrated 
that a significant variation from Zr 02 stoichiometry is not observed; instead a sharp transition 
occurs from Z 1 O 2 to ZrO over a distance of a few nm, with a more gradual change from ZrO to 
oxygen-saturated a-Zr occurring over a distance of tens of nm. 

1.3 Oxide phase stability 

Various investigations have demonstrated the presence of high compressive stresses (0.2-3 GPa) 
in the oxide layer in the plane of the oxide-metal interface [221 - ]?o| . due to the Pilling-Bedworth 
ratio of 1.56 on converting from zirconium to ZrCh- Experimental work had shown that applied 
hydrostatic stress can stabilise the metastable tetragonal ZrC >2 phase due to it having a smaller 
volume than the monoclinic phase (251127) . Recent synchrotron x-ray diffraction work performed 
by Ortner et al. [28j has demonstrated a significant variation in the stress through the oxide layer, 
with higher in-plane compressive stresses observed close to the metal-oxide interface. Similar 
work performed by Polatidis et al. also identified a discontinuity in the stress profile, with a 
sudden increase in compressive stress at a depth within the oxide layer, which had previously 
been identified, using electron microscopy, as the approximate point of transition from monoclinic 
to tetragonal rich Z 1 O 2 [23]. 

Small grain sizes (typically less than 30 nm) can also result in tetragonal phase stability 
with no additional environmental factors |171l29ll30] . Qin et al. [3lj combined these theories 
into a thermodynamic model, suggesting that in any given environment there is a critical grain 
size below which the tetragonal phase is stable. Applied stress increases this critical grain size 
allowing more tetragonal phase to be stabilised. 

Dopant stabilisation can occur due to the incorporation of larger cations that expand the 
oxide lattice, or by the incorporation of lower valence cations, which stabilise the tetragonal 
and cubic fluorite phases by the incorporation of oxygen vacancies. While the concentrations of 
dopant elements required for complete cubic phase stabilisation (e.g. 20 mol% yttrium [32]) are 
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far higher than the concentration of alloying elements in any common zirconium based alloys, 
concentrations as low as 2 mol% yttrium can stabilise the tetragonal phase [33]. When the other 
complementary stabilisation methods (stress and grain size) are present, it is possible that even 
1-2% of trivalent and divalent cations could have a stabilising effect. 

1.4 Tin as an alloying addition 

Tin was added initially to reduce the detrimental effects of nitrogen and carbon impurities, and 
is still present in Zircaloy-2 and Zircaloy-4 as it also improves strength and creep resistance [3J. 
It has, however, been demonstrated that the removal of tin from modern cladding alloys can 
increase the time until transition occurs and lower HPUF [351135] . 

Unlike many other alloying elements, the amount of tin included is lower than the solid 
solubility limit in a-Zr, so that a homogeneous distribution throughout the cladding matrix is 
expected. Wei et al. [35] observed that in the oxide layers of various Zr-Sn-Nb alloys, those with 
lower tin content exhibited a smaller proportion of pre-transition tetragonal phase. A reduction 
in the corrosion rate, due to a delayed transition, was also observed, implying a correspondence 
with tetragonal phase fraction, tin content and corrosion rate, which had not previously been 
suggested. 

In this study, the behaviour of tin as a substitutional and interstitial defect in tetragonal 
Z 1 O 2 is investigated using computer simulation. Quantum mechanical calculations using den¬ 
sity functional theory are performed and comparisons are made to previous experimental data 
concerning the corrosion of tin-containing zirconium alloys. 

2 Methodology 

Simulations were performed using the density functional theory based CASTEP 6.11 code m- 
Ultra-soft pseudo potentials with a cut-off energy of 550 eV were used throughout. The Perdew, 
Burke and Ernzerhof [40] parametrisation of the generalised gradient approximation was em¬ 
ployed to describe the exchange correlation function. A Monkhorst-Pack sampling scheme [U; 
was used for the integration of the Brillouin Zone, with a minimum k-point separation of 
0.045 A -1 . The simulations employed density mixing using the Pulay method W2\ . 
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The energy convergence criterion for self-consistent calculations was set to 1 x 10~ 8 eV and 
the maximum allowed forces between ions was f x fO -2 eV /A. All simulations were performed 
until a maximum difference in energy of 1 x 10~ 5 eV and atomic displacement of 5 x 10 -4 A 
between iterations was achieved. 

Non-defective structures were relaxed under constant pressure to the above convergence crite¬ 
ria. All defective structures were generated from pre-relaxed non-defective structures, and were 
energy minimised under constant volume (cell parameters constrained to maintain the shape 
and volume of the perfect supercell) in order to approximate dilute conditions. A supercell was 
formed from 3x3x2 repetitions of the tetragonal Zr 02 unit cell in the x, y and z directions re¬ 
spectively. This resulted in a 108 atom supercell, which offered a reasonable compromise between 
reducing finite size effects and computation time. 

2.1 Defect formation energies 

The defect formation energies ( E l ) of intrinsic defects were calculated using Equation [1] where 
a point defect X with charge q is formed in a perfect cell. 

E { x q = — ^perfect + Q(Ey BM + He) + -^MP (1) 

i 

fgSct * s the ener gy °f the perfect cell, is the energy of the defective cell, -EVbm is the 

valence-band maximum (VBM) of the perfect supercell, rq is the number of atoms added/removed, 
Hi is the chemical potential of the defect species added/removed. Emp is an energy correction 
calculated using the Makov-Payne method (33] to account for the electrostatic self-interaction of 
defects caused by the use of periodic boundary conditions and a finite supercell size. He is the 
chemical potential of electrons relative to Evbm; taking values of He at the VBM and conduction- 
band minimum, the defect formation energy can be plotted as a function of He across the band 
gap of the material [431135] . 

2.2 Chemical potentials 

For a given set of conditions, the sum of the chemical potentials (/z) per formula unit of the 
constituent species equals the total Gibbs free energy of the solid Z 1 -O 2 , leading to the following 
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relationship, 

MZr0 2 (s) = VZr{T,Po 2 ) + P0 2 ( T ,P0 2 ) ( 2 ) 

Where T is the temperature and po 2 is the oxygen partial pressure. In order to avoid the well- 
known difficulties regarding the inability of DFT to accurately describe the O 2 molecule, we 
adopt an approach first used by Finnis et al. Ml- This avoids the necessity to calculate the 
oxygen chemical potential in DFT by instead using the known experimental formation energy of 
the Z 1 O 2 oxide, and calculating the chemical potential of oxygen at standard temperature and 
pressure using the following relationship: 

A Gf r ° 2 = PZr0 2 (s) — ^Zr(s) ~ P0 2 (T° ,Po 2 ) (3) 

We assume the Gibbs free energy is independent of temperature for the solid species, however for 
the gaseous O 2 , this assumption does not hold true. Instead the oxygen chemical potential at the 
temperature of interest po 2 (T,po 2 ) is extrapolated from the standard state p,o 2 {T°,Pq 2 ) using 
the ideal gas relationship. The full expression for the oxygen chemical potential is as follows: 

Po 2 (T,po 2 ) = po 2 (T°,Po 2 ) + A piT) + \k B log (^j (4) 

and the rigid-dumbell ideal gas for A p(T) can be given by: 

A p(T) = -i (S° 02 - Cp) (T - T°) + CpTlog (5) 

Where Sq 2 is the molecular entropy at standard temperature and pressure and Cp is the constant 
pressure heat capacity of oxygen gas. Values for these two properties were obtained from the 
literature with = 0.0021 eV/K and C£ = 7 k B = 0.000302 eV/K 07]. 

Tetragonal Zr02 is not stable under standard conditions, so the free energy of formation was 
calculated by adding the energy difference between monoclinic and tetragonal Zr 02 (as calculated 
at 0 K by DFT) to the experimentally determined value for the free energy of formation of 
monoclinic Zr02, this resulted in a value of —11.41 eV. 
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2.3 Brouwer Diagram 


The sum of all defects each multiplied by their charge must equal zero, since there is no overall 
charge on the crystal. This can be expressed as follows: 




N c exp 


Eg \ 

k B T ) 


+ N v exp 



= 0 


( 6 ) 


where the first term is the sum of the charges of all ionic defects, the second term is the electron 
concentration and the third term the hole concentration in the crystal. N c and N v are the density 
of states for the conduction and valence bands and E g is the band gap of the crystal. Tetragonal 
Zr 02 is an insulating material and so the concentrations of electrons and holes are expected 
to be sufficiently low that Boltzmann statistics are appropriate. The formation energy for an 
electron in the conduction band will be slightly lower than the value calculated by E g — /j, e due 
to self trapping of electrons. However, in a wide band gap insulator such as ZrC >2 the difference 
in energy will be minimal and so this is an acceptable approximation. 

Using this relationship, the electron chemical potential that ensures charge neutrality in the 
system can be calculated for any given elemental chemical potential. Using this, the concentration 
of individual defects can be calculated and by plotting these as a function of oxygen partial 
pressure, a Brouwer diagram is constructed. 


3 Results and Discussion 

3.1 Interstitials 

Potential interstitial sites were chosen by considering the Wyckoff positions of the space group 
for tetragonal Z 1 O 2 (P^/nmc). The 2b and 4c positions are very similar, however both were 
included to investigate whether there was a preference to be closer to a zirconium cation (as in 
the case of the 4c position) or a position equidistant between two (as in the 2b position). An 
interstitial ion was placed in each non-identical site and an overall charge was applied to the 
supercell corresponding to the charge state of the interstitial ion; -2 to 0 for oxygen, 0 to +4 
for zirconium and tin. A list of the non-identical sites considered, along with the corresponding 
lattice co-ordinates is reported in Tabic [[] 
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A successful convergence was only achieved for oxygen in the 2b location, suggesting that the 
other sites are not stable for any charge state. Convergence was achieved for zirconium and tin 
interstitials on all sites; for both the 2b site was preferred. Tin and zirconium interstitials placed 
at the 4c site migrated to the 2b site during the simulation and the 4d site was not energetically 
favourable for any interstitials, generally exhibiting a formation energy several eV higher than 
the 2b site. 

Table 1: A list of the non-identical Wyckoff positions for tetragonal ZrC >2 considered as possible 
interstitial sites. 


Position 

a 

b 

C 

2b 

0 

0 

0.5 

4c 

0 

0 

-0.25 

4d 

0.5 

0 

0 


3.2 Defect formation energies 

The defect formation energies for single isolated defects as a function of /i e calculated using 
Equation |T] are shown in Figure [T] Oxygen defects exist almost entirely in the fully charged or 
zero charged states, implying that the singly charged state is not thermodynamically favourable 
for either a vacancy or interstitial. 

Zirconium vacancies exist in the -4 oxidation state across the majority of the band gap, with 
a transition through all charge states to 0 observed close to the VBM. The zirconium interstitial 
is also observed in the fully charged +4 state across the majority of the band gap, transitioning 
to +2 and then 0 close to the conduction band minimum (CBM). These observations agree well 
with previous simulation work investigating tetragonal Zr02 I48fl50| . 

Tin interstitial defects are predicted to have a very high formation energy when compared to 
other possible defects and so can be assumed to exhibit a negligible concentration. The formation 
energies for the substitutional Snz r defect are shown in Figure [T| The S n^ r defect is favoured 
across the majority of the band gap, with a transition to S n Zr occurring close to the CBM. The 
S n' Zr defect is not observed at any point in the band gap, this is expected as it would otherwise 
require tin to be in the +3 charge state which would result in an energetically unfavourable 
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Fermi Level n e (eV) Fermi Level p e (eV) 

Figure 1: Formation energies of intrinsic defects (a) and tin substitutional defects (b) as calculated by 
Equation |T| using DFT calculated energies from a reference state of Zr( s ) and plotted as a function of 
the Fermi level (/r e ), from the VBM across the experimental band gap of 5.75 eV 1511 . 

unpaired s-orbital electron. 

3.3 Brouwer diagram 

Following the methodology discussed in Section 12.31 Brouwer diagrams were constructed for 
intrinsic tetragonal Z 1 -O 2 and for tin doped ZrC^, and are shown in Figure [2] As discussed in 
Section 11.31 tetragonal ZrC >2 is stabilised in the undoped oxide layer by a combination of grain 
size and compressive stress. In this work, stress was not applied and so the Brouwer diagrams 
were plotted at 1500 K, the approximate temperature at which tetragonal phase is stable under 
standard conditions, rather than at a normal reactor operation temperature of around 600 K. In 
this regard we follow the approach used in previous DFT studies on the tetragonal system Em¬ 
in the intrinsic case (Figure (2Jr) at very low oxygen partial pressures the uncharged oxygen 
vacancy appears. However, across the majority of the oxygen partial pressures considered the 
dominant structural defect is the fully charged oxygen vacancy, charge compensated by electrons. 
Only at the highest oxygen partial pressures does the fully charged zirconium vacancy, charge 
compensated by holes, begin to appear. 

The diagram that includes substitutional tin (Figure [5 }d) indicates that Sn£. is dominant 
across the majority of oxygen partial pressures considered, with a transition to S n Zr charge 
compensated by Vq' occurring at partial pressures below 10 -31 atm. 
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Figure 2: Brouwer diagram showing the concentrations of point defects in tetragonal Zr02 as a function 
of oxygen partial pressure at 1500K. (a) is the intrinsic system, (b) contains tin at a concentration of 
1 x 1(T 4 at.%. 


3.4 Tin interaction with other alloying elements 

Tin has been shown to have little effect on the corrosion resistance and no measurable effect on 
HPUF in Zr-Sn binary alloys [53] • As previously discussed, the tin Brouwer diagrams (Figure [2]) 
suggest that tin exists as Sn^ T . until extremely low oxygen partial pressures. Given that Sn 4+ 
and Zr 4+ occupy the same lattice site, have similar ionic radii and exist in the same oxidation 
state it is perhaps unsurprising that the Sn^, defect has little measurable effect on the chemistry 
of the system. However, Sn-Zr binary alloys are not used in any reactor applications and so we 
must consider the effect of further alloying additions. 

Niobium containing alloys have been developed with excellent corrosion and HPUF resistance 
and recent experimental work has suggested that the removal of tin, included for its positive effect 
on the mechanical properties of the alloys, improves the corrosion resistance still further [35]. 
Given that tin in isolation has no effect on the corrosion and HPUF properties, it seems plausible 
that there is an interaction between tin and niobium in the alloys. 

Niobium is generally assumed to exist in the N6 5+ state within the oxide layer. If we assume 
the 5+ state, this results in a predominant N b' Zr defect, which agrees with DFT simulations by 
Otgonbaatar et al. [54]. This would lead to the injection of positive charge into the oxide layer. 
Nevertheless, this assumption has been challenged recently through XANES work which instead 
suggests the charge state is between 2+ and 4+ [sang. Recent theoretical and experimental 


10 


















work by Couet et al. [57] has suggested another possible source of positive charge close to the 
metal-oxide interface due to a space charge effect caused by the non-equilibrium distribution of 
electrons in the insulating oxide layer. 

In order to consider the effect of the positive charge, an additional defect was added to the first 
term of Equation [G] This defect was given a charge of +1, and a Brouwer diagram was plotted 
for a range of concentrations. The inclusion of an additional positive defect at a concentration of 
5 x 10 -4 (Figure [3]i) causes the concentration of V^" to dramatically increase and also lowers the 
concentration of Vq* and Vq (the concentration of Vq is unchanged since it is not charged and 
therefore unaffected by charge difference in the system). Importantly, although the concentration 
of S n Zr is increased, S n^ r remains the dominant defect at partial pressures above 10 -31 atm. 

Increasing the concentration of applied charge further (Figure [3]b) causes S n Zr to become the 
dominant tin defect up to a partial pressure of 10~ 20 atm. This increase of around 10 orders 
of magnitude happens at a critical charge concentration and as demonstrated by Figure [4] even 
large changes in charge above or below this critical value have little effect on the crossover point. 
The critical charge at which transition occurs appears to be strongly dependent on temperature, 
with an increase of 500 K causing the value to increase by an order of magnitude, suggesting 
that the Sn^ defect is more favourable at higher temperatures. The existence of a critical value 
has previously been predicted in sapphire by Lagerlof and Grimes [58], who observed that a 
critical doping level caused the concentration of oxygen vacancies to change by 10-20 orders of 
magnitude. Interestingly, as the temperature of the Brouwer diagram calculation was increased 
the change in partial pressure observed at the critical doping level reduced. 
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Figure 3: Brouwer diagrams containing tin at a concentration of 1 x 10 4 at.% and an additional defect 
of charge +1 (q = 1) inserted into the calculation via Equation [6] (a) c = 5 x 10 —4 , (b) c = 2 x 10 -3 . 



0.0001 0.001 o.oi 

Artificial defect concentration (per Zr0 2 ) 


Figure 4: A chart showing the oxygen partial pressure at which the Sn^ r /Sn Zr transition occurs as 
a function of the applied additional charge; 1250 K at a concentration of (2.4 — 2.5) x 10~ 4 , 1500 K at 
(1.15 - 1.18) x 10 -3 , 1750 K at (4.5 - 4.6) x 10 -3 


3.5 Implications for corrosion resistance and HPUF 

It is generally assumed that niobium improves corrosion resistance by suppressing the formation 
of oxygen vacancies in the bulk oxide, thereby reducing the oxygen ion conductivity. This is 
achieved due to the positive charge introduced by the dominant N b* Zr defect, which suppresses 
V** and Vq since they are also positively charged. However, as demonstrated in Figures [3] and |H 
the presence of a concentration of positive charge in the system causes the S n^ r defect to become 
S n Zr . Considering the oxide layer on a Zr-Sn-Nb type alloy, this tin behaviour would negate 
the oxygen vacancy suppression delivered by niobium by charge compensating N b' Zr defects, 
suggesting a possible mechanism by which tin diminishes the corrosion resistance of Zr-Sn-Nb 
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alloys. 

As mentioned in the introduction, oxygen vacancies can have the effect of stabilising the 
tetragonal phase. Thus, a further interesting implication is that the increased Vq* concentration 
in Sn-containing Nb-Zr alloys could result in the increased tetragonal phase fraction observed 
in XRD work by Wei et al. [55] , As corrosion progresses, the metal-oxide interface will move 
further away from the S n Zr substitutional defects interface due to the thickening of the oxide 
layer, which will have the effect of increasing the oxygen partial pressure in the oxide surrounding 
the defects. As demonstrated by Figures [3] and [2 there is a critical partial pressure at which 
the S n Zr defect will reduce to Sra£ r . As this reduction occurs, the increased V“ concentration 
associated with the presence of S n Zr is no longer expected. Thus, the stabilisation mechanism 
causing the increased tetragonal phase volume will no longer exist and a transformation to 
monoclinic phase will occur. The tetragonal to monoclinic phase transformation is associated 
with a volume increase of approximately 4% [59) and so cracking of the oxide layer is expected. 
This process could lead to an earlier onset of transition than would otherwise occur, a conclusion 
that agrees closely with experimental observations that the reduction of Sn content delays the 
first transition significantly G3- 

4 Conclusions 

DFT simulations were used to investigate the defect properties of intrinsic and Sn-doped tetrago¬ 
nal Zr02- The intrinsic case agreed well with previous DFT work performed by Youssef et al. [48] . 
Tin is predicted to exhibit a 4+ charge state and exist almost entirely in uncharged S n^ r , tran¬ 
sitions to a 2+ charge state and thus S n Zr only at oxygen partial pressures below 10~ 31 atm. 

An additional positively charged defect was included in the Brouwer diagram calculations to 
account for the effect of positively charged defect species, in particular N b m Zr . It is predicted that, 
up to a critical concentration, the additional defect has little effect on the dominance of S n^ r 
defect, with the only observed effect being the suppression of oxygen vacancies and increased 
concentration of zirconium vacancies, as would be expected in a Nb-containing system. Above 
a temperature dependent critical concentration, however, the additional positive charged defect 
promotes the transition from S n^ r to S n Zr to occur at an oxygen partial pressure many orders 
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of magnitude higher. The concentration at which this change occurs and the magnitude of the 
change in partial pressure has a strong temperature dependence; a higher charge concentration 
was required at high temperatures but resulting in a smaller shift in partial pressure. 

It is suggested that this change in tin oxidation state will inhibit the oxygen vacancy sup¬ 
pression, which would otherwise be caused by the additional positive charge, thereby accounting 
for the increased corrosion rate and shorter time to transition observed in Zr-Sn-Nb type alloys 
when compared to similar alloys containing niobium but no tin. With no niobium content, we 
expect tin to maintain an overall neutral charge (Sn^. r ) and therefore have little effect on the 
corrosion performance or time to transition. This provides a model framework by which we can 
interpret experimental work on Zn-Sn binary alloys |53] . 
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